{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "73697504",
   "metadata": {},
   "source": [
    "This notebook works as a visualization of the different functions defined in the file `ellipse.py`. It also works as small unit-tests, to check the functionality of these functions. "
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "f34c8632",
   "metadata": {},
   "source": [
    "### Fit ellipse to points"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "86ea4f2f",
   "metadata": {},
   "outputs": [],
   "source": [
    "import ellipse\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7601d8a5",
   "metadata": {},
   "outputs": [],
   "source": [
    "npts = 250\n",
    "tmin, tmax = np.pi / 6, 4 * np.pi / 3\n",
    "x0, y0 = 4, -3.5\n",
    "ap, bp = 7, 3\n",
    "phi = np.pi / 4\n",
    "# Get some points on the ellipse (no need to specify the eccentricity).\n",
    "x, y = ellipse.get_ellipse_pts((x0, y0, ap, bp, phi), npts, tmin, tmax)\n",
    "noise = 0.1\n",
    "x += noise * np.random.normal(size=npts)\n",
    "y += noise * np.random.normal(size=npts)\n",
    "\n",
    "coeffs = ellipse.fit_ellipse(x, y)\n",
    "print('Exact parameters:')\n",
    "print('x0, y0, ap, bp, phi =', x0, y0, ap, bp, phi)\n",
    "print('Fitted parameters:')\n",
    "print('a, b, c, d, e, f =', coeffs)\n",
    "x0, y0, ap, bp, phi = ellipse.cart_to_pol(coeffs)\n",
    "print('x0, y0, ap, bp, phi = ', x0, y0, ap, bp, phi)\n",
    "\n",
    "plt.plot(x, y, 'x')  # given points\n",
    "x, y = ellipse.get_ellipse_pts((x0, y0, ap, bp, phi))\n",
    "plt.plot(x, y)\n",
    "plt.show()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "20723f46",
   "metadata": {},
   "source": [
    "### Project point to ellipse"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "66268a17",
   "metadata": {},
   "source": [
    "This is to project the OCR recognitions to the ellipse"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "54360cdb",
   "metadata": {},
   "outputs": [],
   "source": [
    "x0, y0 = 0, 0\n",
    "ap, bp = 7, 4\n",
    "phi = 0\n",
    "# Get some points on the ellipse (no need to specify the eccentricity).\n",
    "x_e, y_e = ellipse.get_ellipse_pts((x0, y0, ap, bp, phi))\n",
    "plt.plot(x_e, y_e)\n",
    "plt.scatter(x0, y0)\n",
    "\n",
    "# point to project on to the ellipse\n",
    "x = 5\n",
    "y = 10\n",
    "point = np.array([x, y])\n",
    "plt.scatter(x, y)\n",
    "\n",
    "projected_point = ellipse.project_point(point, (x0, y0, ap, bp, phi))\n",
    "\n",
    "x_proj, y_proj = projected_point\n",
    "plt.scatter(x_proj, y_proj)\n",
    "plt.scatter(x_proj, y_proj, marker = 'x')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "7b1b730e",
   "metadata": {},
   "source": [
    "### Intersect line and ellipse"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "faa50404",
   "metadata": {},
   "source": [
    "This is to calculate the intersection of the needle and the ellipse."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b29f4830",
   "metadata": {},
   "outputs": [],
   "source": [
    "line_coeffs = np.array([1,1])\n",
    "line = np.poly1d(line_coeffs)\n",
    "x_start, x_end = 0.1, 2\n",
    "x = np.array([x_start, x_end])\n",
    "\n",
    "x0, y0 = 1, 1\n",
    "ap, bp = 1, 2\n",
    "phi = 0\n",
    "# Get some points on the ellipse (no need to specify the eccentricity).\n",
    "x_e, y_e = ellipse.get_ellipse_pts((x0, y0, ap, bp, phi))\n",
    "plt.plot(x_e, y_e)\n",
    "plt.plot(x, line(x), color='orange')\n",
    "\n",
    "intersection_point = ellipse.get_line_ellipse_point(line_coeffs, x, (x0, y0, ap, bp, phi))\n",
    "print(intersection_point)\n",
    "x = intersection_point[0]\n",
    "y = intersection_point[1]\n",
    "\n",
    "plt.scatter(x, y, marker = 'o', c='red', s=100)\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "1050c260",
   "metadata": {},
   "source": [
    "### Get middle between two points on ellipse"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "17f5aab6",
   "metadata": {},
   "source": [
    "This is to calculate the point between the start and end point."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "56cfb63f",
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(5,5))\n",
    "\n",
    "start_point = np.array((-0.5,-2))\n",
    "end_point = np.array((1.5,-1))\n",
    "\n",
    "plt.scatter(start_point[0], start_point[1], marker = 'o', c='blue', s=100)\n",
    "plt.scatter(end_point[0], end_point[1], marker = 'o', c='green', s=100)\n",
    "\n",
    "x0, y0 = 1, 1\n",
    "plt.scatter(x0, y0, marker = 'o', c='black', s=100)\n",
    "ap, bp = 1, 2\n",
    "phi = np.pi*5/4\n",
    "ellipse_params = (x0, y0, ap, bp, phi)\n",
    "# Get some points on the ellipse (no need to specify the eccentricity).\n",
    "x_e, y_e = ellipse.get_ellipse_pts(ellipse_params)\n",
    "plt.plot(x_e, y_e)\n",
    "\n",
    "theta_start = ellipse.get_polar_angle(start_point, ellipse_params)\n",
    "theta_end = ellipse.get_polar_angle(end_point, ellipse_params)\n",
    "\n",
    "start_proj = ellipse.get_point_from_angle(theta_start, ellipse_params)\n",
    "end_proj = ellipse.get_point_from_angle(theta_end, ellipse_params)\n",
    "plt.scatter(start_proj[0], start_proj[1], marker = 'o', c='blue', s=100)\n",
    "plt.scatter(end_proj[0], end_proj[1], marker = 'o', c='green', s=100)\n",
    "\n",
    "zero_point = ellipse.get_point_from_angle(0, ellipse_params)\n",
    "plt.scatter(zero_point[0], zero_point[1], marker = 'o', c='red', s=100)\n",
    "\n",
    "theta_middle = ellipse.get_theta_middle(theta_start, theta_end)\n",
    "middle = ellipse.get_point_from_angle(theta_middle, ellipse_params)\n",
    "\n",
    "plt.scatter(middle[0], middle[1], marker = 'o', c='violet', s=100)\n",
    "    \n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "dbba4910",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "ffd351a3",
   "metadata": {},
   "source": [
    "# Warp ellipse"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "b01932c7",
   "metadata": {},
   "source": [
    "Here is some code to warp the ellipse, such that it becomes a circle"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "01b299c5",
   "metadata": {},
   "outputs": [],
   "source": [
    "import cv2\n",
    "import numpy as np\n",
    "\n",
    "def warp_ellipse_to_circle(image, ellipse_center, ellipse_axes, ellipse_angle):\n",
    "    image_height, image_width = image.shape[:2]\n",
    "\n",
    "    # Define the source points (the coordinates of the four corners of the ellipse)\n",
    "    x, y = ellipse_center\n",
    "    major_axis, minor_axis = ellipse_axes\n",
    "    source_points = np.array([[-major_axis/2, -minor_axis/2],\n",
    "                            [major_axis/2, -minor_axis/2],\n",
    "                            [major_axis/2, minor_axis/2],\n",
    "                            [-major_axis/2, minor_axis/2]], dtype=np.float32)\n",
    "    square_size = max(major_axis, minor_axis)\n",
    "    destination_points = np.array([[-square_size/2, -square_size/2],\n",
    "                            [square_size/2, -square_size/2],\n",
    "                            [square_size/2, square_size/2],\n",
    "                            [-square_size/2, square_size/2]], dtype=np.float32)\n",
    "\n",
    "\n",
    "    rotation_matrix = np.array([[np.cos(phi), np.sin(phi)], [-np.sin(phi), np.cos(phi)]])\n",
    "    source_points = source_points @ rotation_matrix\n",
    "    destination_points = destination_points @ rotation_matrix\n",
    "\n",
    "    source_points[:, 0] += x\n",
    "    source_points[:, 1] += y\n",
    "    destination_points[:, 0] += x\n",
    "    destination_points[:, 1] += y\n",
    "\n",
    "    source_points = source_points.astype(np.float32)\n",
    "    destination_points = destination_points.astype(np.float32)\n",
    "\n",
    "    # Calculate the perspective transformation matrix\n",
    "    transformation_matrix = cv2.getPerspectiveTransform(source_points, destination_points)\n",
    "\n",
    "    # Warp the image\n",
    "    offsetSize=0\n",
    "\n",
    "    warped_image = cv2.warpPerspective(image, transformation_matrix, (image_width+offsetSize, image_height+offsetSize))\n",
    "\n",
    "    return warped_image, transformation_matrix\n",
    "\n",
    "# Example usage\n",
    "image_path = '/home/$USER/$image_path'  # Replace with the path to your image\n",
    "image = cv2.imread(image_path)\n",
    "\n",
    "ellipse_center = (355, 208)  # Replace with the coordinates of the ellipse center\n",
    "ellipse_axes = (200, 80)  # Replace with the major and minor axes lengths of the ellipse\n",
    "\n",
    "x0 = ellipse_center[0]\n",
    "y0 = ellipse_center[1]\n",
    "ap = ellipse_axes[0]\n",
    "bp = ellipse_axes[1]\n",
    "phi = np.pi/4\n",
    "x, y = ellipse.get_ellipse_pts((x0, y0, ap, bp, phi))\n",
    "plt.figure()\n",
    "plt.imshow(image)\n",
    "plt.plot(x, y)\n",
    "\n",
    "plt.show()\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "15171584",
   "metadata": {},
   "outputs": [],
   "source": [
    "warped_image, transformation_matrix = warp_ellipse_to_circle(image, ellipse_center, ellipse_axes, phi)\n",
    "plt.figure()\n",
    "point_warp = [440,125]\n",
    "plt.imshow(warped_image)\n",
    "plt.scatter(point_warp[0], point_warp[1])"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "9fdb0828",
   "metadata": {},
   "source": [
    "This is code to remap a point in the warped image to the original image"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "68809e7d",
   "metadata": {},
   "outputs": [],
   "source": [
    "def map_point_original_image(point_warp, transformation_matrix):\n",
    "    inverse_transformation_matrix = np.linalg.inv(transformation_matrix)\n",
    "\n",
    "    point_in_original_image = cv2.perspectiveTransform(np.array([[point_warp]], dtype=np.float32),\n",
    "                                                    inverse_transformation_matrix)[0][0]\n",
    "    return point_in_original_image\n",
    "\n",
    "point_in_original_image = map_point_original_image(point_warp, transformation_matrix)\n",
    "\n",
    "plt.figure()\n",
    "plt.imshow(image)\n",
    "plt.scatter(point_in_original_image[0], point_in_original_image[1])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "09733c3a",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "855a6c2a",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "image_height, image_width = image.shape[:2]\n",
    "\n",
    "# Define the source points (the coordinates of the four corners of the ellipse)\n",
    "x, y = ellipse_center\n",
    "major_axis, minor_axis = ellipse_axes\n",
    "source_points = np.array([[-major_axis/2, -minor_axis/2],\n",
    "                            [major_axis/2, -minor_axis/2],\n",
    "                            [major_axis/2, minor_axis/2],\n",
    "                            [-major_axis/2, minor_axis/2]], dtype=np.float32)\n",
    "\n",
    "rotation_matrix = np.array([[np.cos(phi), np.sin(phi)], [-np.sin(phi), np.cos(phi)]])\n",
    "source_points = source_points @ rotation_matrix\n",
    "\n",
    "source_points[:, 0] += x\n",
    "source_points[:, 1] += y\n",
    "\n",
    "ellipse_points_x, ellipse_points_y = ellipse.get_ellipse_pts((x0, y0, ap, bp, phi))\n",
    "plt.figure()\n",
    "plt.imshow(image)\n",
    "plt.plot(ellipse_points_x, ellipse_points_y)\n",
    "plt.scatter(source_points[:,0], source_points[:,1])\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b9c7b2d9",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the destination points (the coordinates of the four corners of a square)\n",
    "square_size = max(major_axis, minor_axis)\n",
    "destination_points = np.array([[-square_size/2, -square_size/2],\n",
    "                            [square_size/2, -square_size/2],\n",
    "                            [square_size/2, square_size/2],\n",
    "                            [-square_size/2, square_size/2]], dtype=np.float32)\n",
    "\n",
    "\n",
    "\n",
    "destination_points = destination_points @ rotation_matrix\n",
    "destination_points[:, 0] += x\n",
    "destination_points[:, 1] += y\n",
    "\n",
    "plt.figure()\n",
    "plt.imshow(image)\n",
    "plt.plot(ellipse_points_x, ellipse_points_y)\n",
    "plt.scatter(destination_points[:,0], destination_points[:,1])\n",
    "\n",
    "source_points = source_points.astype(np.float32)\n",
    "destination_points = destination_points.astype(np.float32)\n",
    "\n",
    "\n",
    "# Calculate the perspective transformation matrix\n",
    "transformation_matrix = cv2.getPerspectiveTransform(source_points, destination_points)\n",
    "\n",
    "# Warp the image\n",
    "warped_image = cv2.warpPerspective(image, transformation_matrix, (image_width, image_height))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eac197f2",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "01ac2f61",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9dfffdfd",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
